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ABSTRACT 

The effect of photon frequency redistribution by line branching on mass-loss rates for 
hot luminous stars is investigated. Monte Carlo simulations are carried out for a range 
of OB-star models which show that previous mass-loss calculations which neglect non- 
resonance line scattering overestimate mass-loss rates for luminous O stars by ^ 20 
per cent. For luminous B stars the effect is somewhat larger, typically ~ 50 per cent. 
A Wolf-Rayet star model is used to investigate line branching in the strong wind limit. 
In this case the effect of line branching is much greater, giving mass-loss rates that 
are smaller by a factor ^ 3 from computations which neglect branching. 
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1 INTRODUCTION 

It has been known for over thirty years that luminous, 
hot stars have powerful stellar winds (Morton 1967). These 
winds are important in a variety of astrophysical contexts, 
most notably because of their dynamical influence on the in- 
terstellar medium, their effect on stellar evolution and their 
role in metal enrichment of galaxies. 

The proposal that these winds are driven by radiation 
pressure, in particular by photon scattering in spectral lines, 
was developed by Lucy & Solomon (1970) and is now gen- 
erally accepted (Kudritzki & Puis 2000). Since the iden- 
tification of the mass-loss mechanism, attention has been 
concentrated on the development of techniques which can 
predict reliable mass-loss rates. Currently, the most accu- 
rate mass-loss rates for hot stars are obtained from Monte 
Carlo simulations of radiation transport in winds (de Koter, 
Heap & Hubeny 1997; Vink, de Koter & Lamers 1999, 2000; 
hereafter VdKL99, VdKLOO respectively), an approach to 
mass-loss calculations first applied by Abbott & Lucy (1985, 
hereafter AL85). In this method, the global constraint of 
energy conservation is used to determine the mass-loss rate 
algebraically while the wind terminal velocity is fixed at the 
observed value. This is in contrast to other techniques, such 
as 'CAK' theory (Castor, Abbott & Klein 1975), where the 
mass-loss rate is an eigenvalue determined at a critical point 
and the terminal velocity is calculated. Monte Carlo simula- 
tions have also been used to study mass-loss as a function of 
metallicity (Vink, de Koter & Lamers 2001) and mass-loss 
from luminous blue variables (Vink & de Koter 2002). 
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A weakness of previous determinations of mass-loss 
rates from Monte Carlo simulations is the assumption that 
all line scattering events may be treated as resonance scat- 
terings. This assumption was made by AL85 and has been 
followed by others in subsequent work (including the mass- 
loss calculations of de Koter et al. 1997 and VdKL99,00). 
It has the advantage of significantly simplifying the compu- 
tation and can be justified for many of the strong lines, 
even some non-resonance lines, in several important ions 
(see AL85). However, although remarkably successful for 
computing mass-loss rates, the assumption is not valid in 
general: Lucy (1999b) included downwards branching in a 
Monte Carlo spectral synthesis code for supernovae and 
found this to have a significant impact on the emergent spec- 
trum. 

Following recent developments of the Monte Carlo 
energy-packet technique (Lucy 2002, 2003), it has become 
feasible to lift the resonance scattering assumption and al- 
low branching to occur in line scattering events. Branch- 
ing may effect mass-loss rates since it allows photon leakage 
from one part of the spectrum to another (see e.g. the dis- 
cussion by Owocki & Gayley 1999). If only resonance scat- 
tering is permitted, the frequency of a given photon packet 
is gradually redshifted in the frame co-moving with the stel- 
lar wind (which accelerates outwards in the star's reference 
frame): photons have no opportunity to jump in frequency. 
Branching alters this by allowing photons to be re-emitted 
in a different part of the spectrum following a line scatter- 
ing event. The net effect of line branching on many photons 
is to preferentially redistribute photons from regions of the 
spectrum with many strong spectral lines (e.g. ultraviolet) 
to relatively sparse regions (e.g. infrared). This will reduce 
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the net line force since photons undergo fewer line scatter- 
ings. The purpose of this paper is to quantify this effect and 
determine the extent by which the assumption of resonance 
scattering causes mass-loss rates to be overestimated. It is 
to be expected that the magnitude of the effect will depend 
on the typical number of scatterings that photons undergo 
as they traverse the wind: in the limit where no photons 
are scattered more than once (single scattering limit) the 
frequency of re-emission is irrelevant to the mass-loss calcu- 
lation meaning that line branching can have no effect. 

Most of the results presented in this paper relate to 
OB stars, the stars for which VdKLOO have computed mass- 
loss rates. However, calculations appropriate for Wolf-Rayot 
(W-R) stars are also presented. The winds of W-R stars 
present the greatest challenges to the theory of radiation 
driven mass-loss, requiring significantly higher efficiency in 
the transfer of momentum from the radiation field to the 
matter than OB star winds. Lucy & Abbott (1993, here- 
after LA93) investigated line driving in W-R stars and con- 
cluded that it was possible to obtain wind performance num- 
bers T] = ^Vooc/L„ which are significantly in excess of one, 
as required for W-R winds ($ is the mass-loss rate, Voo is 
the terminal velocity of the wind and L, is the luminos- 
ity of the wind- free star), de Koter et al. (1997) performed 
Monte Carlo simulations with relatively sophisticated mod- 
els for specific WN5h and Of/WN stars (as classified by 
Crowther & Dessart 1998) to compute line-driven mass-loss 
rates for comparison with values determined from ffa obser- 
vations. They found that the computed mass-loss rates were 
too small by factors <~ 2 (note, however, that there are 
significant uncertainties in the mass-loss rates determined 
from observations) . A complete understanding of the mass- 
loss mechanism in W-R stars does not yet exist, but it seems 
certain that line and continuum driving have some role (see 
e.g. Nugis & Lamers 2002; Lamers & Nugis 2002). This in- 
vestigation is of particular relevance to W-R stars because 
the high density and ionisation stratification of W-R winds 
means that multiple scattering is much more prevalent than 
in OB staxs. Since the effect of line branching on the radia- 
tive acceleration can only be important when multiple scat- 
tering occurs, W-R stars should exhibit the most extreme 
effects. 

Subsequent to the commencement of the work presented 
here, Onifer & Gayley (2003) have published a complemen- 
tary study of frequency redistribution in stellar winds. They 
present a paramctcriscd analytical theory appropriate for 
an idealised spectrum consisting of two frequency domains, 
one line-rich and the other line-poor. However, they did not 
perform any calculations with a real line list meaning they 
were unable to obtain quantitative estimates of the effect 
of frequency redistribution on mass-loss rates. In contrast, 
the numerical computations that are presented here do not 
require any idealisations of the distribution of line opacity 
across the spectrum and use a real line list, allowing quan- 
titative results to be derived. 

Study of line branching is also relevant for predictions 
of ionising fluxes from hot stars. Calculations which neglect 
line branching, for example CoStar models (Leitherer et al. 
1996, Schaerer & dc Koter 1997), may overestimate the ion- 
ising flux as a result (Crowther et al. 1999). The new Monte 
Carlo techniques employed here (Lucy 2002,2003) will read- 
ily allow this effect to be quantified in the near future. 



In Section 2, the stellar wind models that will be used 
for the mass-loss calculations are discussed. The Monte 
Carlo approach used for computing mass-loss rates is de- 
scribed in Section 3, and Section 4 introduces the atomic 
data used in the calculations. The results for models of OB 
stars are presented in Section 5 and for a W-R star model 
in Section 6. Conclusions are drawn in Section 7. 



2 MODELS 

2.1 Velocity and density 

As discussed above the objective of this paper is a differential 
study of the effects of line branching and not the develop- 
ment of high quality stellar wind models. Therefore, simple 
"core-halo" models are adopted for the sake of tractability. 
They have been constructed following AL85 and LA93. 
The one-dimensional velocity law 

y{r)=Vo + {v^-Vc)[l-^y (1) 

is adopted where Rc is the core radius (the base of the wind) 
and Vc is the wind velocity at Re- The adopted values of Rc, 
Vc, Voo and l3 are given in Sections 5 and 6. The mass density 
p is determined from the equation of mass conservation 

4iTr^pv = $ . (2) 

2.2 Ionisation and excitation 

The approach to ionisation and excitation recommended by 
Lucy (2003) is adopted here. Namely, only approximate non- 
LTE solutions are pursued and analytic formulae are used 
for both ionisation and excitation. This approach has been 
used in many previous studies (including AL85 and, in part, 
VdKL99,00) and is particularly suitable for this study since 
it decouples the Monte Carlo simulation of the radiation 
field from the computation of level populations. This pro- 
vides two advantages: first the need for iteration when con- 
structing models is removed and, more importantly, it allows 
direct study of the effect of line branching on the radiation 
field in isolation, without any feedback influence on the level 
populations. 

For ionisation, the modified nebular approximation is 

adopted 




where Nj is the number density of the jih ionisation stage, 
Ne is the free electron number density, W is the dilution 
factor which is defined in terms of the photospheric radius 
(Rp) following LA93, i^j is the fraction of recombinations 
to ion j that go direct to the ground state, Te is the elec- 
tron temperature and Tr is the radiation temperature. The 
last bracket in equation (3) is the ionisation fraction com- 
puted using the Saha equation at temperature Tr. The par- 
tition functions used for this calculation are non-LTE, hav- 
ing been computed using diluted Boltzmann population ra^ 
tios for non-metastable levels (see below). The C- values used 
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by AL85 have been adopted throughout this paper. The ra- 
diation and electron temperatures used for the calculations 
are discussed in Sections 5 and 6. 

When determining level populations, following AL85 
and VdKL99, a distinction is made between mctastablc lev- 
els and those with permitted electric dipole decay routes. For 
a metastable level i, the population relative to the ground 
state 1 is assumed to be given by the Boltzmann formula 

^=(^) (4) 

while for all other levels it is assumed that the population 
is driven by radiation such that 

!^ = w('-^] . (5) 

Equations (3), (4) and (5) allow all level populations to 
be determined for a given value of p. Although approximate, 
Springmann & Puis (1998) have demonstrated that in a typ- 
ical O star, these formulae provide excellent approximations 
to ion/level populations obtained from much more detailed 
non-LTE calculations. 



3 THE RADIATION FIELD 

The radiation field is simulated by a Monte Carlo calcu- 
lation using the Macro Atom approach developed by Lucy 
(2002). The principles of the method are discussed by Lucy 
(2002,2003). The details of the Monte Carlo simulation used 
here are presented below, with particular emphasis on the 
differences between the calculations discussed by Lucy and 
those performed here. 

3.1 Discretization of model 

For the Monte Carlo calculation, the wind models (Section 

2) are discretized into a series of thin spherical shells. The 
velocity is allowed to vary in accordance with equation (1) 
throughout the model, but all other parameters are taken to 
be constant within each shell. To evaluate the density (p), for 
each shell equation (2) is used for twenty different values of r 
which are equally spaced in a:: = Rc/r between the inner and 
outer boundaries of the shell. The adopted density for the 
shell is the a;"''- weighted mean of these twenty values. (This 
weighting is chosen to preserve mass in spherical geometry.) 
The dilution factor (W) is evaluated at the mid-point of 
each shell. The innermost shell lies on the surface of the 
core r = Rc, the outermost shell extends to r = 100-Rc- 
Typically 100 shells of equal thickness in x = Rc/r are used 
in the calculations. 

3.2 Creation of r-packets 

The Monte Carlo quanta used in the simulation are indivis- 
ible packets of radiative energy, termed r-packets by Lucy 
(2002). The Monte Carlo method imposes radiative equilib- 
rium by conserving both the number of packets and their 
co-moving energy throughout the calculation. Thus packets 



are not created in the domain of the simulation but are only 
launched into the inner shell (from the surface of the star) 
and propagate until they either re-enter the star through the 
lower boundary of the innermost shell or emerge from the 
upper boundary of the outermost shell. 

The radiation field incident from the stellar core is mod- 
elled by a black-body spectrum at the temperature of the 
stellar surface Te{Rc). A lower limit to the wavelength of 
the incident spectrum is imposed at Amin = 228 A on the 
assumption that higher energy photons are immediately re- 
moved via bound- free absorption by He ii. 

The frequencies i> of the incident r-packets are equally 
spaced in Inz^, and their energies are assigned according to 
the Planck spectrum. Packets are created with wavelengths 
up to lOOAmin- The initial direction cosines ^ are chosen 
randomly, assuming no limb-darkening. Typically, 2 x 10^ 
packets are used in the calculation. 



3.3 Propagation of packets 

As the packets propagate through the wind they undergo in- 
teractions with the material in the wind. Immediately after 
their creation at the lower boundary, after a packet crosses a 
shell boundary and also after any interaction with wind ma- 
terial, the optical depth that a packet will traverse before its 
next interaction is chosen at random using t = — In « where 
is a random number in the interval to 1. Contributions 
to the optical depth along the packet's flight path from both 
continuous and line opacity arc then considered to determine 
whether the packet undergoes a continuum scattering, a line 
absorption or whether the packet will cross a boundary into 
another shell. The method for selecting which event occurs 
follows closely that described by Lucy (2003), sec. 5, and the 
reader is referred to that paper for a full explanation. The 
subsections that follow briefly explain how the continuum 
and line contributions to the opacity are computed. 



3.3.1 Electron scattering 

The only continuum process included in the calculations is 
scattering by free electrons. Although bound-free processes 
are important in spectral synthesis, they are not the domi- 
nant contribution to the radiation driving force in OB stars 
and are therefore not important in this differential study of 
line branching. 

Electron scattering is assumed to be isotropic and co- 
herent in the co-moving frame. 



3.3.2 Line absorption 

The large velocity gradients in a stellar wind mean that line 

processes can be treated in the Sobolov limit. Lucy (2002) 
gives expressions for the Sobolev radiative rates and optical 
depth for the case of an atmosphere in homologous expan- 
sion (his equations 20, 21 and 22). However, for the case 
of non-homologous expansion in a stellar wind, these ex- 
pressions must be generalised. The radiative rates for spon- 
taneous emission and absorption in the transition between 
upper level i and lower level j are given by the standard 
expressions 
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Rij = AijPijUi and Rji = {BjiUj — Bijni)Jj 



(6) 



where A and B are the usual Einstein coefHcients, the n's 
are level populations (number densities) and Jji is the in- 
tegrated mean intensity in the line. The escape probability, 
Pji is given by (Klein & Castor 1978) 



Pi 



2i-i 



and the Sobolev optical depth Ts {ji) is 



d/i 
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4n dv/ds 



where Vij is the line frequency and 



dv 
ds 



c 



[1 - H )- + H — 
r dr 



(7) 



(8) 



(9) 



is the rate of change of co-moving frequency with path- 
length for a photon with rest frame frequency ur at ra- 
dius r travelling in the direction specified by /x. The de- 
pendence of Ts on di^/ds docs not significantly complicate 
the Monte Carlo calculation since di//ds may be taken as 
constant along any given photon trajectory within a shell, 
provided that the shells are sufficiently thin. 

In this formulation, the interaction between a packet 
and a spectral line with which is comes into resonance is 
described by a single scattering event. In reality, when a 
photon enters resonance with a line it may scatter many 
times in that line before it is redshifted out of the region in 
which interaction occurs (the resonance zone). Since ther- 
malisation is not included here, these multiple interactions 
with a particular line are readily described by a single scat- 
tering event and the usual Sobolev escape probability (see 
de Koter et al. 1997 for a discussion of thermalisation during 
multiple interactions within a resonance zone) . Accordingly, 
throughout this paper the term single scattering is used to 
refer to photon packets that interact with one line during 
the Monte Carlo simulation, even though this interaction 
may physically consist of several scatterings. Multiple scat- 
tering refers to energy packets that undergo two or more 
independent interactions with different spectral lines during 
their fiight through the model. 

The Sobolev optical depth of a line is used to determine 
whether or not a packet is absorbed when it comes into res- 
onance with the line. If absorption occurs the Macro Atom 
method developed by Lucy (2002) is used to determine the 
atomic level from which re-emission occurs and a random 
number is used to determine in which of the transitions from 
the level the packet is emitted, following Lucy (2003). The 
Macro Atom formalism does allow for coUisional, bound-free 
and free- free processes but, for simplicity, these are not in- 
cluded here. 

The calculation of the Macro Atom upward jumping 
probabilities requires J to be known for each spectral line. 
Lucy (1999b, 2003) has shown how efficient Monte Carlo es- 
timators may be constructed for integrals of the radiation 
field. Following Lucy, the estimator used here for a transi- 



tion between upper level i and lower level j in a particular 
shell p is 



Jji — 



AirAt Vr 
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(10) 



where At is the time interval represented by the Monte Carlo 
simulation and Vp is the volume of the shell p. The sum runs 
over all packets which come into resonance with the line 
in shell p during the simulation, e is the packet energy in 
the co-moving frame and Ts is the Sobolev optical depth of 
the transition. Note that since Ts is angle dependent every 
packet encounters a different optical depth: thus the fac- 
tors involving Ts must be placed inside the summation. For 
the branching calculations presented in the later sections 
of this paper, the Jji-values are calculated iteratively for 
each model, with initial values taken as a dilute black-body 
radiation field. The values of Jji were found to converge to 
sufficient accuracy for the mass-loss calculations within only 
one or two iterations. 

For all the models discussed in this paper, in addition 
to finding mass-loss rates when branching is permitted, cal- 
culations have been performed in which the Macro Atom 
formalism has not been employed, and all lines have been 
treated as resonance lines (no line branching). The compar- 
ison of results using the full Macro Atom approach and the 
resonance scattering assumption forms the basis of the dif- 
ferential analyses presented in Sections 5 and 6. 

Once the emission frequency of the packet is deter- 
mined, the new direction of propagation is determined by 
randomly sampling the distribution function for the escape 
probability. Thus the probability of emission with direction 
cosine n is 

1 _ p-Ts(M) 

P(m) dM = K dM (11) 

TsW 

where K is a normalisation constant. 



3.4 The mass-loss rate 

Once all the packets have propagated through the model and 
have either escaped from the upper boundary or re-entered 

the stellar core, the rate of energy deposition in the wind 
can be computed. It is given by AL85 



Lt — 



J2k 



J2iU 



(12) 



where the e's are the rest frame packet energies and the sums 
run over all packets as they enter the calculation (index k), 

escape the wind (index i) and re-enter the star (index j). 
This energy goes into lifting material out of the gravitational 
field of the star and so the mass-loss rate implied by the 
Monte Carlo calculation is given by (LA93) 



2Lt 



MC 



(13) 



where the escape speed at r = Rc- This can be com- 

pared to the mass- loss rate assumed in the modelling ($). 
By performing calculations for models with different values 



Table 1. The elements and ions included in the line list. 



Element 


Ions 


Element 


Ions 


H 


I 


He 


I,II 


C 


I - IV 


N 


I - VI 


O 


I - VI 


F 


I - VI 


Ne 


I - VI 


Na 


I - VI 


Mg 


I - VI 


Al 


I - VI 


Si 


I - VI 


P 


I - VI 


S 


I - VI 


CI 


I - V 


Ar 


I - V 


K 


I - V 


Ca 


I - VI 


Ti 


I - VI 


Cr 


I - VI 


Mn 


I - VI 


Fe 


I - VI 


Co 


I - VI 


Ni 


I - VI 







of A = $MC can be found as a function of $ and so the 
self-consistent value of $ (i.e. when A = 1) can be identified. 

This approach docs not look for local dynamical con- 
sistency between the driving force and the velocity law: the 
mass-loss rate is determined algebraically from the global 
requirement A = 1 rather than as an eigenvalue at a critical 
point in the flow. LA93 used this method and showed that, 
for their W-R model, the Monte Carlo radiative line accel- 
eration differed locally from that implied by the assumed 
velocity law by up to 0.29 dex. Ultimately models with re- 
alistic calculations of the line acceleration which impose lo- 
cal dynamic consistency are desirable but, as discussed by 
LA93, in view of the other simplifications made, to pursue 
local consistency here is not warranted. 
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Table 2. Stellar parameters used for ^ Puppis. From Puis et al. 
(1996). 



Parameter 


Value 


Mass, M* 


52.5 M0 


Luminosity, L, 


106 Lq 


Effective temperature, T^g 


42,000 K 


Terminal velocity, Voo 


2250 km s"! 



Boltzmann equation. The electron temperature, Te is as- 
sumed to be constant throughout each model and is set by 

Te = 0.9Teff . (14) 

For their model of Puppis, AL85 adopted a constant 
radiation temperature somewhat lower than the effective 
temperature of the star. A similar approach is used here. 
In principle the radiation temperature for each ion may be 
determined iterativcly during the Monte Carlo simulation 
(see e.g. Lucy 1999a). Iterative calculations of Tr were per- 
formed using the Pup model discussed below in which the 
radiation temperature for each ion was computed using the 
simulated radiation field above its ground state ionisation 
edge. The results of these calculations could be adequately 
described by a constant value of Tr for ionisation edges be- 
low that of He ii and no radiation above this edge owing 
to the truncation of the incident radiation field at 228 A. 
Therefore, to avoid the need for iteration, in the OB star 
calculations that follow 



4 ATOMIC DATA 

Atomic models and radiative transition probabilities are 
needed for the computation of the mass-loss as described 
above. For the calculations presented here, the atomic data 
has all been taken from the Kurucz & Bell (1995) atomic 
database. All elements with solar abundances X/H < 3 x 
10~* have been neglected. For the elements which are in- 
cluded all ionisation stages up to vi that arc in the database 
are used (higher ionisation species are not important for the 
OB-star models below). The elements and ionisation stages 
from which lines have been included in the calculations pre- 
sented in Section 5 are listed in Table 1. 

For each ion a cut-off in log gf was made to exclude 
weak transitions (in order that the line list is of a manage- 
able size). Atomic models were then constructed for each 
ion, based on its line list and a complete, frequency ordered 
line list assembled for the computations. In total, around 
2.5 X 10'' transitions are included in the line list. 



5 OB STARS 

In this section, calculations of mass-loss rates for O and B 
stars will be presented. The models used here have all been 
constructed using a velocity law (equation 1) with /3 = 1, 
Vc = and the density determined from equation (2). In 
all cases the core radius (i?c) is determined from the lumi- 
nosity (L) and effective temperature (Tes) using the Stefan- 



Tr = T^g - 9, 000 K (15) 

is adopted for ious with ionisation potentials less that He ii 
and no ionisation to higher ions is included. The decrement 
of 9, 000 K is larger than the decrement adopted by AL85 

and was chosen in order to give reasonable agreement be- 
tween the models presented here and the more sophisticated 
models used by VdKLOO. 

5.1 C Pup 

The first model that will be used to show the differential 
effect of line branching on mass-loss calculations is chosen 
to represent the well studied 04 I (f) star, ( Puppis. The 
adopted stellar parameters are taken from Puis et al. (1996) 
and listed in Table 2. There is significant uncertainty in 
the mass: the Puis et al. value was derived spectroscopi- 
cally rather than from evolutionary calculations (evolution- 
ary masses are typically higher - e.g. 70 Mq as adopted 
by Lamers & Nugis 2003). The same element composition 
used by AL85 is assumed: i.e. enhanced helium, but all other 
abundances assumed to be solar. 

To investigate the effect of line branching on the mass- 
loss rate of ^ Puppis a range of models with different in- 
put values of 3> have been used to compute A($) (see Sec- 
tion 3.4), initially including the complete treatment of line 
branching. Separate calculations were then performed to 
find A($) when line branching is neglected and all line scat- 
terings are treated as resonant scattering. Fig. 1 shows A($) 
versus $ for both sets of calculations. 
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Figure 1. A as a function of <I> for Puppis: the boxes (connected 
by solid line) show computations with line branching and the 
triangles (connected with dashed line) show computations with 
resonance line scattering only. The dotted horizontal line is A = 
1.0 (self-consistent model) and the vertical dotted line indicated 
the observed mass- loss rate (from Puis et al. 1996). The A values 
are accurate to around ±0.01 dex. 



The calculations yield a self-consistent (A = 1) mass- 
loss rate for Puppis of5.16xl0~^ Mq yr~^ if line branching 
is included and 6.41 x 10^^ M0 yr^^ if only resonance scat- 
tering is permitted. A sufficient number of packets were used 
that these values of 3> are correct to better than ±2 per cent. 
Thus the differential effect of line branching in the C Puppis 
model is a reduction in the mass-loss rate by around 20 per 
cent. 

This reduction occurs as a result of photon leakage: the 
net effect of line branching is the transfer of packets from fre- 
quency ranges with many strong lines to regions with fewer 
lines. This reduces the average number of lines with which 
a packet interacts and thus reduces the mass-loss rate. This 
effect is illustrated by Fig. 2 which shows the number of 
line scatterings that packets undergo when line branching 
is allowed compared to when it is not. Clearly, line branch- 
ing significantly reduces the number of packets that undergo 
multiple (> 10) line scattering events. Of course the wind 
is sufficiently diffuse that, in both cases, the radiative accel- 
eration is dominated by the large number of packets which 
undergo relatively few scattering and so the differential ef- 
fect is relatively small (the wind performance number for 
the C Puppis model is only 77 = 0.56). Clearly, however, the 
effect will be bigger when multiply scattered packets become 
more dominant: i.e. when rj > 1. 

Fig. 3 shows the fractional contribution to the energy 
deposition rate from all the ions included in the (" Puppis 
modelling. These were computed by summing the contribu- 
tions to the energy deposition rate from line scattering for 
each ion during the Monte Carlo simulation, not by per- 
forming separate computations with lines of only one ion at 
a time. The results shown arc from a calculation in which 
branching was included. For ^ Puppis the global mass-loss 
rate has significant contributions from Fe and from several 
of the light elements (S, Ar, N, O, C). When branching is 
not permitted the relative contributions of the elements are 
very similar except that iron contributes very slightly more: 
this is expected since the iron group elements have the most 
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Figure 2. The number of packets versus the number of line scat- 
terings for Puppis from Monte Carlo simulations with a total 
of 2 X 10® packets. The solid histogram shows a simulation where 
line branching is allowed. The dotted histogram is the case where 
only resonance scattering is included. 
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Figure 3. The relative importance of the different elements/ions 
to the mass-loss-rate for ^ Puppis. The chart shows the fraction 
of energy deposited in the wind by each element in a Monte Carlo 
simulation including branching. The fractions axe not significantly 
different if branching is neglected. The bars for each element are 
subdivided by ionisation stage. Note that the fractions do not sum 
to 1.0 because of the contribution of scattering by free electrons 
(not shown). 

complex spectra with many bunches of lines at similar wave- 
lengths meaning that they are the most sensitive to photon 

leakage. 

These results must be put into context regarding their 
relevance to our understanding of massive stars. Firstly, let 
us compare the derived mass-loss rates for 1^ Puppis with 
those derived from models made by others. AL85 adopted 
a slightly smaller stellar mass (40 M0) and higher termi- 
nal velocity (2600 km s^^) than used here and obtained a 
mass-loss rate for ^ Puppis of 5.2 x 10~® M0 yr~^ from 
their standard model, with values ranging from about 4 - 
7 xlO~® Mq yr~^ from their models with different treat- 
ments of ionisation and assumed velocity law. This is in good 
agreement with the mass-loss rate derived here with reso- 
nant scattering only (AL85 did not attempt to include line 
branching in their calculations) . Using the parameters from 
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Table 2, the VdKLOO mass-loss recipe gives a mass-loss rate 
of 9.4 X 10^® M0 yr~^, similar to the value quoted for Pup- 
pis by VdKLOO. This is a little higher than the mass-loss rate 
obtained here, but it is important to note (see Fig. 1) that A 
is a weak function of $ and thus relatively small differences 
in the modelling approach (computation of A) lead to larger 
differences in the derived globally self-consistent value of 
If # = 9.4 X 10~® M0 yr~^ is used in the modelUng here, 
a Monte Carlo simulation (with resonance scattering only) 
leads to A = 0.8, indicating that the models used here re- 
produce the VdKLOO energy deposition rates adequately, to 
within 20 per cent. Given that many simplifying assump- 
tions have been made here (e.g. the stellar radiation field 
is assumed to be black-body, a single value of Tr has been 
adopted for all ions at all radii, no bound-free or free-free 
processes are included, a core-halo approach is adopted in 
the modelling), some discrepancy with the VdKLOO results 
is expected. 

From observations of Ha, Puis et al. (1996) have de- 
rived a mass-loss rate of 5.9 x 10~® yr~^ for ^ Puppis. 
This value is indicated by the vertical line in Fig. 1, and is 
in good agreement with the computed mass-loss rates dis- 
cussed above. However, there is significant uncertainty in 
the observational mass-loss rate: estimates from radio fluxes 
give only 2.4 x 10"*^ M© yr"^ (Lamers & Leitherer 1993). 
Clumping in the wind may affect the accuracy of these deter- 
minations although this effect should be small for OB stars 
(Kudritzki & Puis 2000). In view of the large observational 
uncertainty in the mass-loss rate, the effect of line branch- 
ing can be regarded as small and, for the present, is not 
important for modelling mass-loss from f Puppis to observ- 
able accuracy. However, the effect is not utterly negligible 
and upon a reduction in the observational uncertainty by a 
factor of a few, it should become detectable (note, however, 
that since stellar winds show intrinsic variability at the level 
of several tens per cent, time variability would also need to 
be taken into account in reliable modelling of mass- loss rates 
to such high accuracy). 

5.2 Grid of models 

Following on from the detailed discussion of C Puppis above, 
in this section the effects of line branching on the com- 
puted mass-loss rates are investigated for a range of OB-star 
models spanning roughly the same parameter space as the 
VdKLOO calculations. For each model (defined by the stel- 
lar parameters L»,M*,reff and Voo), two A($) curves have 
been computed, one with line branching included and one 
without. The models considered are listed in Table 3 and 
the A(<1>) curves are shown for each model in Fig. 4. For 
each model, the globally self-consistent (A($) = 1) mass- 
loss rate was found both with and without line branching 
included (#bran and 3>res, respectively). These are given in 
Table 3. The Monte Carlo errors in the mass-loss rates are 
less than 2 per cent. The table also gives the wind perfor- 
mance number (r;) calculated with and without line branch- 
ing, the VdKLOO mass-loss rate (3>vdKL) for comparison, the 
ratio of the branching to resonance-only mass-loss rates and 
the characteristic density of the model p — "l'bran/(47ri?ct'tx)). 
The models are labelled A - O, being ordered first by de- 
creasing effective temperature and secondly by decreasing 



Consider first the high temperature models appropriate 
for O stars (A - H). For these models the A($) curves (Fig. 
4) are all approximately straight in log-log space around 
A = 1. The effect of line branching (i.e. the vertical offset 
between the two curves in Fig. 4) is relatively small for all 
these models - it leads to a reduction in the self-consistent 
mass-loss rate by between three and 25 per cent, relative to 
calculations which neglect branching. At a given temperar 
ture (e.g. compare models C - F) the reduction is closely 
correlated with the characteristic dciisitj' (p). This is ex- 
pected since the effect of line branching is related to multiple 
scattering which is more common at high densities. 

Models I - N all have T^a = 20,000 K (early B stars). 
The first four of these models (I - L) have terminal velocities 
which are typical for stars cooler than the bi-stability jump 
{voo = l-SVgsc where vl^c is the effective escape velocity: 
see VdKL99,00) while the other two have higher terminal 
velocities, characteristic of stars hotter than the jump {voo ~ 
2.6uo5:). For these models, the A(<I>) curves around A = 1 
are more complex than for the higher temperature models 
discussed above (see Fig. 4) : all the low temperature models 
have A = 1 close to a region in which A($) passes through 
a plateau. These plateaus are the result of changes in the 
driving mechanism as a function of the density (which is 
determined by the value of $ adopted in the model) and are 
thus physically related to the bi-stability jumps discussed by 
VdKL99,00. Specifically, the plateau seen here is associated 
with the balance between Fe ill and Fe ll: in models on 
the low-density (i.e. low-$) side of the plateau, Fe iii is the 
single most important ion for driving the wind while on the 
high-density side Fe ll provides the largest contribution. This 
transition in the driving mechanism leads to the complex 
shape of the A($) curve. 

These plateaus influence the effect of line branching on 
the derived mass-loss rates. Consider Models I and J: for 
these models, the vertical offset between A($) computed 
with and without branching is comparable (see Fig. 4). For 
Model I, A = 1 occurs at values of $ beyond the plateau, in 
a region where A($) varies rapidly with $. This leads to a 
relatively small difference between $res and $bran. However, 
for Model ,], A = 1 lies in the region where A($) varies 
very slowly. This means that the difference between $rcs and 
$bran is much greater for this model. Models K, L and M are 
all similax to, but less extreme than, Model J. Conversely, 
Model N behaves similarly to Model I since it has A = 1 
significantly below the plateau. 

Taken together. Models I - N (and also the lower tem- 
perature Model O) indicate that the effect of line branching 
is generally greater at lower temperatures (B stars) than at 
high temperatures (O stars). They also show that subtleties 
in the modelling (in particular changes in the shape of A(3>) 
resulting from changes in the driving mechanism around, for 
example, the bi-stability jumps) can influence the differen- 
tial effect of line branching such that it cannot be described 
by a simple scaling with p and Teff alone. 

For most of the models, the mass-loss rate computed 
with resonance scattering agrees with that found by VdKLOO 
to within a factor of two or better, which is adequate given 
the differential nature of the calculations presented here. 
Those models where the agreement is slightly poorer (Mod- 
els L, K, N and M) all have self-consistent mass-loss rates 
($res) which lie close to plateaus in A($) - since A($) varies 




Figure 4. A($) for models A to O. In each case the boxes (connected with solid lines) show the results with branching and the triangles 
(dashed lines) show the results without branching. The horizontal dotted line shows A = 1.0 (self-consistent model) and the vertical 
dotted line indicated the mass-loss rate from the VdKLOO recipe. 
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Table 3. Computed mass-loss rates for a selection of OB-star models. ^>res and $bran are the mass-loss rates computed with resonance 
scattering only and with line branching, respectively. $vdKL is the mass-loss rate from the VdKLOO recipe. The terminal velocity (uoo) is 
taken as 1.3 times the effective escape speed on the cool side of the bi-stability jump and 2.6 times on the hot side (VdKLOO). 

?7 is the wind performance number, 77 = ^Vooc/L^,. p = ^'bran/C^T-R^ Uoo) is the characteristic mass density in the model. The results are 
accurate to ~ ±2 per cent. 



Model 


log-L, 






«^oo(«oo/fefc) 


log *bran 




log *res 


^res 


log *VdKL 


*bran/*res 


logp 










(km s~~"'") 


(Mr, vr"!) 








(Mr. vr"!) 




( ff cm ~^ ^ 


A 


6.0 


60 


50,000 


2606 (2.6) 


-5.01 


1.24 


-4.89 


1.64 


-5.05 


0.75 


-12.67 


B 


6.0 


80 


50,000 


3265 (2.6) 


-5.23 


0.94 


-5.14 


1.12 


-5.21 


0.82 


-12.99 


C 


6.0 


60 


40.000 


2085 (2.6) 


-5.29 


0.52 


-5.17 


0.68 


-5.04 


0.76 


-13.24 


D 


6.0 


80 


40,000 


2612 (2.6) 


-5.51 


0.39 


-5.42 


0.49 


-5.20 


0.81 


-13.56 


E 


5.5 


40 


40,000 


2660 (2.6) 


-6.24 


0.23 


-6.19 


0.26 


-5.90 


0.89 


-13.79 


F 


5.0 


20 


40,000 


2619 (2.6) 


-6.95 


0.14 


-6.93 


0.15 


-6.61 


0.94 


-14.00 


G 


6.0 


60 


30,000 


1564 (2.6) 


-5.75 


0.13 


-5.67 


0.16 


-5.32 


0.83 


-14.07 


H 


5.0 


20 


30,000 


1964 (2.6) 


-7.03 


0.09 


-7.02 


0.09 


-6.89 


0.97 


-14.45 


I 


6.0 


60 


20,000 


521 (1.3) 


-4.54 


0.72 


-4.37 


1.07 


-4.58 


0.68 


-13.09 


J 


6.0 


80 


20,000 


653 (1.3) 


-5.19 


0.20 


-4.67 


0.68 


-4.75 


0.30 


-13.84 


K 


5.5 


40 


20,000 


665 (1.3) 


-6.27 


0.06 


-5.90 


0.13 


-5.45 


0.43 


-14.43 


L 


5.0 


20 


20,000 


655 (1.3) 


-6.98 


0.03 


-6.77 


0.05 


-6.15 


0.62 


-14.62 


M 


6.0 


60 


20,000 


1042 (2.6) 


-5.78 


0.08 


-5.44 


0.18 


-5.06 


0.46 


-14.63 


N 


5.0 


20 


20,000 


1309 (2.6) 


-7.30 


0.03 


-7.23 


0.04 


-6.63 


0.85 


-15.25 





6.0 


60 


15,000 


391 (1.3) 


-4.57 


0.51 


-4.47 


0.64 


-4.30" 


0.80 


-13.49 



" This value was read from Fig 2. of VdKLOO since the mass-loss recipe (section 5 of VdKLOO) does not fit this model well. 



slowly in these regions relatively minor inadequacies in the 
models (i.e. the computation of A) lead to larger errors in 

the self-consistent value of In any case, these errors are 
of little concern since it can bo seen from Fig. 4 that the 
differential effect of branching (the offset between the pairs 
of curves in the figure) is not a strong function of 3>. 

In conclusion, for hot ~ 40, 000 K stars with lumi- 
nosity ~ 10''Lq calculations which neglect branching over- 
estimate the mass-loss rate by ~ 20 per cent while for 
lower luminosity stars with similar temperatures the effect 
is smaller. The overestimation factor can be larger for the 
cooler (~ 20, 000 K) B stars with similar luminosities (typi- 
cally about 50 per cent and in extreme cases up to a factor of 
a few) and can be significantly affected by variations in the 
ionisation balance (and therefore wind driving mechanism) 
in the models. 



6 WOLF-RAYET STAR 

W-R stars have larger mass-loss rates (up to ~ 
10"'' M0 yr"'^) than OB-stars but comparable luminosities. 
LA93 and Gayley, Owocki & Cranmer (1995) have modelled 
line driving in W-R winds and shown that the high density 
means that multiple scattering is very important in these 
winds (the wind performance number is typically ~ 7, far 
in excess of the single scattering limit). Thus, since the in- 
fiuence of line branching is linked to the predominance of 
multiple scattering, W-R star winds are likely to show a 
large differential effect if line branching is introduced to the 
calculation. This effect is quantified here. 

The W-R wind model adopted for this study is that of 
LA93 (M, = 12.6 Mq, Rc = 2.5 R©, L, = 2.8x lO'" Lq, Vc = 
10 km s~', Vao = 2500 km s~^, /3 = 1). These parameters 
were chosen by LA93 to represent a generic WN5 star. The 
model is similar to the OB-star models discussed above, the 



main differences being the temperature stratification and 
the element abundances. These are discussed briefly below. 

As discussed by LA93 and subsequently confirmed by 
Schulte-Ladbeck, Eenens & Davis (1995) and Herald et al. 
(2000), observations suggest that there is a significant out- 
ward variation of ionisation in W-R winds. Ionisation strat- 
ification has also been predicted theoretically from W-R 
type models (e.g. de Koter et al. 1997). To account for this, 
LA93 imposed a simple variation of radiation temperature 
with velocity (their equation 11) which is also used here. 
Their standard value for the radiation temperature at Doo 
of logTfl = 4.35 is adopted here. The outward variation of 
electron temperature is less critical to the calculation: it is 
obtained from the modified Mihie-Eddington temperature 
distribution (Lucy 1976), following LA93. In addition, the 
usual dilution factor (W) is replaced by the modified dilu- 
tion factor (yVi) defined by LA93. 

The same elemental composition used by LA93 is 
adopted. Thus X = 0.0, Y = 0.98 and Z = 0.02. The car- 
bon, nitrogen and oxygen are taken from a WNE model 
(Maeder & Meynet 1987) and the other metal abundances 
are scaled from solar values to give Z = 0.02. 

As in the OB-star models, the radiation field at the 
base of the model r = _Rc is assumed to be black-body with 
a high energy cut-off at the He 11 ionisation edge. This cut-off 
is inconsistent with the adopted ionisation structure (i.e. if 
there is really no radiation shortward of the He 11 edge then 
there arc no photons to produce the high ionisation stages 
by photoionisation). However, this is of little consequence 
here since in reality there will be a locally produced radia- 
tion field that can produce the high ions. In addition, the 
calculations presented here are relatively insensitive to the 
wavelength cut-off in the incident radiation field since most 
of the input energy is at longer wavelengths: this has been 
checked explicitly with Monte Carlo simulations in which the 
high energy cut-off was placed at higher energy (by a factor 
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Figure 5. A as a function of <I> for the W-R model: the boxes 
(conneeted by solid line) show computations with line branching 
and the triangles (connected with dashed line) show computations 
with resonance line scattering only. The dotted horizontal line is 
A = 1.0 (self-consistent model). The A values are accurate to 
around ±0.01 dex. 



of ten) - the results of those simulations are not significantly 
diflFerent from those presented below. 

The black-body temperature adopted for the incoming 
radiation field is Tcff,c, the effective tornpcraturc determined 
from L, at Re- To allow for the higher temperatures deep in 
W-R winds, the line list used for the OB staxs (see Section 4) 
was expanded to include lines from higher ionisation stages 
(vil and vill) of the most important iron group elements 
(iron and nickel). In the calculations that follow, however, 
these ion stages provide a negligible contribution (<~ 1 per 
cent) to the radiative driving because of the assumed shape 
of the input stellar radiation field - most of the energy is 
at wavelengths that arc too long for the packets to interact 
with the strong lines of these ions. 

Fig. 5 shows A versus $ for W-R models computed 
with and without line branching in the Monte Carlo sim- 
ulation. As expected, Fig. 5 shows that line branching has 
a far larger effect on the mass-loss calculation for a W-R 
star than for OB stars: the difference between A computed 
with and without branching for a given value of $ is much 
greater than for any of the models discussed in Section 5. 
The self-consistent W-R mass-loss rate determined without 
branching is 1.8 x 10~^ Mq yr~^ (which is close to the value 
of 2.1 X 10"^ Mq yr"^ found by LA93). When line branching 
is included the mass-loss rate is only 5.4 x 10~^ Mq yT~^, a 
reduction by more than a factor of 3. Note that while some 
of the cool OB models discussed in Section 5 also predicted 
large reduction factors ('I'rcs/'l'bran), these were not because 
branching had such a large effect on the computed values 
of A for a given value of $ but were because of the relative 
insensitivity of A to the value of $ adopted in the modelling. 

The magnitude of this effect can be explained following 
the discussion of the O stars above. Fig. 6 shows the distri- 
bution of the number of line scatterings for the W-R model. 
The reduction in the number of packets that are scattered 
many times as a result of line branching is dramatic in the 
W-R model: without branching, individual packets undergo 
up to ~ 700 line scatterings (in a simulation with 10® pack- 
ets) but with branching included no packets undergo more 
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Figure 6. The number of packets versus the number of line scat- 
terings for the W-R model from Monte Carlo simulations with 
a total of 10*' packets. The solid histogram shows a simulation 
where line branching is allowed. The dotted histogram shows the 
case where only resonance scattering is permitted. 

than 90 line scattering events. Comparing Figs. 2 and 6, 
shows that multiple scattering is much more important in 
the W-R model than in the ( Puppis model. Thus photon 
leakage from one part of the spectrum to another has a much 
greater influence on the derived mass-loss rate. 

Fig. 7 shows the contribution of the different ions to the 
mass-loss rate for the W-R model. The lower panel is from 
a simulation with resonant scattering only and the upper 
shows calculations with line branching. The W-R mass-loss 
rate is dominated by the iron group elements (Fe and Ni). 
When branching is included the relative contribution of the 
iron group elements drops and the lighter elements become 
more important: as discussed in Section 5.1, this is expected 
because of the greater complexity of the iron group ions. 
The effect here is more dramatic than in ^ Puppis, but is 
still small (the fractional contribution of Fe and Ni to the 
Une driving is smaller when branching is allowed by only 
~ 10 per cent). 

There is some observational evidence to suggest that 
the wind acceleration is more gradual in W-R stars than 
OB stars: e.g. based on infrared data Ignace, Quigley & 
Cassinelli (2003) have recently suggested /3 values of 2-3 for 
a WN6 star, and they do not rule out larger values. There- 
fore, to investigate the sensitivity of the results presented 
here to the adopted value of (3, calculations have been per- 
formed using the same W-R parameters as above except 
that (3 = 3 was used for the velocity law. In this case, the 
derived mass-loss rates were smaller: 7.2 x 10~^ Mq yr~^ 
and 4.4 x 10^'' Mq;, with and without line branching, re- 
spectively. The differential effect of line branching is smaller 
here than for the (3 = 1 case, but is still significant (a factor 
of ~ 1.6). 



7 CONCLUSIONS 

In the previous sections, simple wind models for a range 
of hot stars have been used to investigate the effect of line 
branching on mass-loss rates. For O stars, where accurate 
mass-loss rates are important for evolutionary calculations. 
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Figure 7. The relative importance of the difFcreiit cleineiits/ioiis 
to the mass-loss-ratc for the W-R modeh The upper panel shows 
the contributions when branching is included and the lower panel 
shows the case where only resonance scattering is included. In 
both cases the chart shows the fraction of energy deposited in the 
wind by each element. The bars for each element are subdivided 
by ionisation stage. Note that although the seventh and eighth 
stages of iron and nickel are included in the simulation they make 
a negligible contribution which is not shown here. The fractions 
do not sum to 1.0 because of the contribution of scattering by 
free electrons (not shown). 



the overestimation is relatively small, typically < 25 per 

cent and is correlated with the predominance of multiple 
scattering. For these stars the difference in the mass-loss 
rate computed with and without line branching is compa- 
rable to or less than typical observational errors in 0-star 
mass-loss rates. Thus, despite the assumption of resonant 
scattering, in most cases the VdKLOO mass-loss calculations 
are not significantly in error. In the future, however, when 
observational constraints become tighter it will be necessary 
to include line branching in realistic mass-loss calculations. 

The effect of line branching in relatively cool (~ 
20, 000 K) luminous stars (early B stars) is significantly 
larger than in the hotter stars (on average about 50 per 
cent), but it is also more sensitive to the details of the mod- 
elling. 

Line branching is most important under conditions typ- 
ical of a W-R star. The W-R calculations presented here 
suggest that neglecting line branching results in an overes- 
timation of the mass- loss rate by a factor ~ 3 for /? = 1. 
The effect is less severe if a larger value of /3 is adopted but 



is still significant at least for /3 ~ 3. This has implications 
for modelling line driving in W-R stars: branching makes it 
significantly harder for line driving to accelerate the winds 
of these stars. There are many unanswered questions relat- 
ing to W-R stars, particularly regarding the roles of lino and 
continuum driving throughout their winds. But the results 
presented here clearly show that branching significantly in- 
fiuences mass-loss calculations which involve multiple scat- 
tering and that future studies of line driving in W-R stars 
must incorporate these effects. 
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